function Y = CCD_peak(X)
% Small function to detect cosmic ray peaks in CCD images (hot pixels)
% NOTE: It doesn't work very well, better if we use some spatial filtering
% Syntax:
% Y = CCD_peak(X)
% Input:
% X = Initial Image 
% Output: 
% Y = Corrected Image
% EGC 2009/04/30
n = 1:numel(X);
[nX mX] = size(X);
X = reshape(X,[1 numel(X)]);
% Y = X;
peak = fpeak(n,X,30,[n(1),n(end),max(X)/2,inf]);
X = reshape(X,[nX mX]);
% for k=1:length(peak(:,1)),
%     Y(peak(k,1)) = mean (X([peak(k,1)-100:peak(k,1)-1 peak(k,1)+1:peak(k,1)+100]));
% end
% [zmax,imax,zmin,imin] = extrema2(X);

% [j k] = ind2sub(size(X),peak(:,1));
[j k] = ind2sub([nX mX],peak(:,1));

% pix = [34.3245614035088 42.0906432748538;
%     34.3245614035088 47.7748538011696;
%     40.1959064327485 47.7748538011696;
%     40.219298245614 41.9035087719299];
% Y = roifill(X,pix(:,1),pix(:,2));
q = 3;
for n = 1:length(j),
    c = [k(n)-q; k(n)-q; k(n)+q; k(n)+q];
    r = [j(n)-q; j(n)+q; j(n)+q; j(n)-q];
    Y = roifill(X,c,r);
end
% pixel(44,37) & (45,37)